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1 Introduction 


Spectral methods are characterized by the representation of the solution to a differential 
equation in terms of a truncated series of smooth global functions which are known as trial or 
basis functions. The basis functions are usually chosen to be the eigenfunctions of a singular 
Sturm-Liouville problem (Gottlieb and Orszag, 1977). It is this choice which is responsible 
for the superior approximation properties of spectral methods over other standard methods 
of discretization. For linear problems possessing smooth solutions these eigenfunctions yield 
expansions that converge asymptotically faster than any finite power of AM. 

Two areas of research in spectral methods which are receiving much attention at the 
current time are the construction and analysis of well-posed approximations to the Stokes 
and Navier-Stokes equations and the development of methods which can be applied easily 
to problems defined in complex domains. With respect to the first, it is well-known that 
in the primitive variable formulation the velocity and pressure approximation spaces need 
to be compatible to avoid problems of ill-conditioning. This is similar to the Babuska- 
Brezzi condition required for the corresponding finite element approximation spaces. In two 
dimensions it is possible to avoid this difficulty by reformulating the governing equations in 
terms of a stream function. The governing equation is then fourth-order, and nonlinear in 
the case of the Navier-Stokes equations. In this paper we seek to construct pseudospectral 
approximations to fourth-order differential equations with the ultimate goal of applying them 
to solve the nonlinear stream function formulation of the Navier-Stokes equations. 

Secondly, the development of techniques for handling complex geometries is essential 
if spectral methods are to be applied to problems defined in more than just the simplest 
domains. The basic idea behind domain decomposition is to break up the domain into smaller 
simpler subdomains in which spectral approximations can be used. The approximations are 
suitably linked by appropriate interface continuity conditions. The way in which this is 
implemented is important if the full power of the spectral method, in terms of the accuracy 
of the approximation, is to be achieved. 

In this paper we shall restrict ourselves to the model fourth-order problem in one and 
two dimensions. Starting from a variational formulation of the problem we shall derive a 
corresponding collocation problem complete with interface conditions. In a domain decom- 
position setting this approximation will be chosen to be C 1 continuous implicitly. In addition 
C continuity across the subdomain boundaries is achieved asymptotically as the order of 
the approximation is increased. 

Although there are many applications of spectral methods to solve second-order elliptic 
partial differential equations in the literature there is little previous work on fourth-order 
problems even though the regularity of the solution to these problems is generally higher than 
for second-order problems. Some interesting ideas are proposed in the works of Morchoisne 
(1984) and Orszag (1971). Bernardi and Maday (1988) give a survey of strategies that may 
be employed for fourth-order problems. 

Maday and Metivet (1986) have studied Chebyshev spectral and pseudospectral approx- 
imations of the stream function formulation of the Navier-Stokes equations. They prove the 
convergence of the schemes and derive error estimates in weighted Sobolev spaces Kara- 
georghis and Phillips (1989a, 1989b) use a spectral collocation strategy to solve for the laminar 
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flow through a channel contraction again using a stream function formulation for moderate 
values of the Reynolds number. They use a domain decomposition method to subdivide 
the flow region into rectangular subdomains and patching to piece the solutions together, m 

some sense, across the subdomain interfaces. 

In a collocation method the choice of the collocation points is crucial. In spectral methods 
they are always chosen to be the nodes of a Gaussian quadrature rule principally for two 
reasons. First, the Lagrange interpolating polynomial which interpolates data at these nodes 
has good approximation properties. Secondly, the collocation method may be shown to be 
equivalent to a variational formulation of the problem when the same Gaussian quadrature 
rule is used to approximate the integrals appearing in this formulation. For second-order 
problems the Gauss- Lobatto nodes are used because the boundary conditions can then be 
imposed efficiently. This leads to an optimal error in the resulting spectral approximation 
(Canuto et al (1987)) For fourth-order problems two boundary conditions are imposed on 
the solution. These are usually of Dirichlet and Neumann type. The imposition of these 
boundary conditions is facilitated by the construction of a generalized Lagrange interpolating 
polynomial which interpolates the function at the interior nodes and the function and its 
derivative at the boundary nodes. The generalized Gaussian quadrature rule associated wi 
this interpolating polynomial can then be derived. Quadrature rules of this form are quite 
well-known in the theory of numerical integration (see, for example, Golub and Kautsky 
(1983), and the references therein). Golub and Kautsky (1983) describe how the weights 
in these quadrature rules may be determined computationally. In this paper closed form 
expressions for the weights are derived using the properties of orthogonal polynomials 

We show that, for fourth-order problems, the natural choice of nodes are the zeros of cer- 
tain Gegenbauer (or ultraspherical) polynomials. Explicit representations for the quadrature 
weights are derived for evaluating integrals of the form 

f w\(x)f(x)dx , 


where the weight function takes the form 


,2\A 


w\(x) = (1 - x 2 ) 


A> -1. 


The particular form of these weights is given when A = 0 (the Legendre weight function) and 
\ = —1/2 (the Chebyshev weight function). The interior nodes in the case when — / 

are the zeros of T%{x) whereas the interior Gauss-Chebyshev-Lobatto nodes are the zeros of 

TN Allocation scheme for solving a fourth-order model problem is derived by considering 
a variational formulation of the boundary value problem with suitably defined inner prod- 
ucts. The two formulations are shown to be equivalent if the inner product in the discrete 
variational problem is defined by the generalized Gauss quadrature rule. The linear system 
of equations which derives from this collocation scheme is ill-conditioned. The condition 
number of the coefficient matrix scales like 0{N S ) where N is the order of the approxi- 
mation. An efficient preconditioner for this system based on a low order finite difference 
approximation to the same differential operator is presented. The combination of general- 
ized Gaussian quadrature rules with spectral methods has also been propose y ernar i 
et al. (1990). This idea is extended to multi-domain problems in the present paper. rseu- 
dospectral approximations which are C x continuous at the subdomain interfaces are used to 
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approximate the solution in each subdomain. The discrete variational problem enables us 
to derive interface continuity conditions which, in the asymptotic limit, result in C 3 con- 
tinuous approximations. The variational formulation is used to provide an analysis of the 
collocation scheme for domain decomposition. The analysis shows that the pseudospectral 
approximation is optimal in the sense that it is of the same order as the corresponding error 
in the best approximation. Numerical results are presented in which the usual exponential 
convergence behaviour of spectral approximations is exhibited. Finally, the extension of 
the method to two dimensions is described and numerical results presented for a number of 
model problems. An application of the method to the solution of the biharmonic equation 
in a non-rectangular domain, an L-shaped region, for which standard spectral methods are 
not applicable is presented. 


2 Variational Formulation of the Model Problem 


In this section we consider the variational formulation of the fourth-order model problem. 
Consider the fourth-order boundary-value problem 


d 4 u 

IP = /(l) ’ 


-1 < X < 1 , 


u(±I) = 0, 


£( ±1 ) = °' 


(1) 


where f(x ) is a given source function. It is well-known that, for any / <E tf -2 (-l,l), (1) 
has a unique solution u G ( see Grisvard (1985), for example ). A collocation 

scheme for solving (1) is derived by considering a variational formulation of the problem 
with suitably defined inner products. 

To set up the variational formulation we need to define function spaces for each A > — 1 . 
Let L 2 (— 1, 1) be the Hilbert space defined by 




{ v : (—1, 1) — > R is measurable 
fl j w\(x)v 2 -(x)dx < oo 


endowed with the inner product 


) 


(u, v)\ = j ^ w\(x)u(x)v(x)dx. (2) 

We also introduce the Sobolev space H 2 (— 1, 1) defined by 

Hl {~ 1, !) = < v : ^ € L a(-1, 1) , 0 < k < 2} 

with corresponding norm 


Let HIq(- 1, 1) be the subspace of 1) defined by 

// a %(-1 > 1 ) = {u€// a 2 (-1,1): v(±l) = 0; v'(±l) = 0 }. 

Consider now the bilinear form a A (.,.) defined on 

//’(-!,!) x^aV-M) 

by ,i 

a A (u, v ) = J x 

For any / € 1) the fourth-order model problem (1) is equivalent to the variational 

problem : Find u € f/ Ai o( — 0 suc b that 

a,\(u,v) = (/, v) A , Vc e 1). ( 4 ) 

Bernardi and Maday(1991) have proved the following result. 

Proposition 2.1 Let A satisfy -1 < A < 1. The bilinear form a A is elliptic on , 1) 

and continuous on // A ( — 1,1) x Hf 0 ( — 1,1), Lt. 

a II Vue//J„(-l.l), (5) 

| a\(u,v) |< ft || u lla.xll v || 2 ,a , V u £ H\{— 1> 1) > v € V, 0 (—l, 0 ’ 

respectively , where a and ft are positive constants. 

This is a generalization of an earlier result of Maday (1990) for A = -1/2. An immediate 
consequence of Proposition 3.1 is the following theorem: 

Theorem 2.1 Let A satisfy -1 < A < 1. For any f € H?(- 1,1) the variational problem 
(4) has a unique solution u € H\ Q { — 1,1)- Moreover, it satisfies 

II U || w < C II f ||h- j (-i,i) • ( ? ) 

3 Pseudospectral Approximation 

We consider the pseudospectral discretization of the fourth-order problem (1). Let Pn(- 1, 1) 
denote the space of algebraic polynomials of degree N or less on the interval [-1, 1]. Letjr;, 
l<j<N-\beN- \ distinct points in the interval [-1, 1] with x x = -1 and xjv-i - 1. 
Throughout this paper we take N > 4 in order to have at least one point in the interior of 
the interval [—1, 1]. Suppose that the values fj of some function f(x) are given at the points 
x - , l < j < N — 1, together with the values f{ and f' N _ x , of f'{x) at x = Xj, and x = xyv-i, 
respectively. To set up the pseudospectral approximation of (1) which automatically satisfies 
the boundary conditions it is necessary to construct the Lagrange polynomials for this data. 
Define the polynomials u(x) and II(x) by 

v(x) = (1 - x 2 ) 2 , TT(x)=n(x-x J ) (8) 

J=2 
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g‘ V e a n n by e Verifled tha ‘ ‘ he U ® range P0,3 ' n0mial » f d <*“* N which interpolates this data is 

AT-l 

Pn{x) — ^ fjhj(x) -f /jAi(a:) + /Ar_jftw-i(:r), 


( 9 ) 


where 


fij(x) = i 


_n’(r|v(i) 


(x-x ; )f T (*,)„(*,) > 2 < J < N - 2 , 

y = i,Ar-i, 


and 


ft (x) = JJM. ffdfl 
"(*i) 3(*;) ’ 


».f T \_ n (*H*) v(x) 

A [ -(«-»,)» • 

The corresponding integration rule based on these points 

Lx w *( x )f( x ) dx = S Wj/fo) + ^/'(x,) + WtV-l/'^JV-O , 

i= i 


(10) 

( 11 ) 

( 12 ) 

(13) 


can be shown to be exact for all / € P w _ 3 (-1, 1) if the interior nodes Xj ,2<j<N~2 
are c osen o e the zeros of the Gegenbauer or ultraspherical polynomial G<)+ 2 > of degree 

when we investigate the solution of differential equations we will only consider the c*se A=1l 

P s'rr' S ^ Me giVe " in the A PP“dix 

L _ i 2 n u! ! “u j f j" te USed t0 denote equation (m) from the Appendix 
( fhi'p’olynil 1 !: yfl d <T< ™ i“ MdrT"! be “7? V «.owing PP 

polynomials in /u'rn we ^ «* 


w 


’ ~ Lx w ^ x ) h i( x ) dx , l<j<N-l, 

’ j = Lx w ^ x %( x ) dx , j = 1 , N ~ 1 . 


(14) 

(15) 

re,eVa r the «*"*"* P * P « h. that they 
of the’reprJ^nXTLr o wf htf “^7 P ” Me,n detaiU 

un^own'^dL valin' th7IZ metnC SyS ‘ em ° f liDear *° >* solved fofthe 

(9 imar„bt^ deriVe ( r ° ri ® inal reSUlt in Which ex P' idt ^presentations for the weights 
(9,10) are obtamed us.ng the properties of the ultraspherical polynomials (see the Applndix) 
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Let us begin with the weights u>,-, 2 < J < N - 2, associated with the interior points. The 
polynomials TT(rr) and <#!?(*) are related by 


g&2(*> 

n(x) - 


(16) 


since they are of the same degree and have the same zeros, where A N - a is the leading 
coefficient of G*^>(x). Thus using (10) and (14) we may wnte 


W; = 


’i ~ r (W< 

^N - 3 1 


1 Z 1 WJ 


«a(«)»(*)giE?(») 


(x - Xj) 


dx, 2<j<N-2. 


(17) 


In order to determine the value of this integral, we make use of the Christoffel-Darboux 
identity: 


n -3 < 7 ( , A+2) (x')Gi A+2) (y) G^?(aO Gffi*?(y) - G N-l\ x ) G N- 2 ] (y}^i_ 
— - — - — — = * T TZ I7\ Aw_ 

^Vl 


7n- 3(* — y) 


-3 

As--i 


(18) 


* (A+2) 

where 7t , (1 < * < N -3) is defined by (A.4). Now we replace y by z, where G %tl (*i) = 0 
then (18) reduces to 

ggg(f) _ _ or’woriiai. do) 

7JV-3-/4jV_2 ( x ~ x i) fc=o 

If we now multiply both sides of (19) by (1 - x 2 )*+ 2 <^( x ) and integrate with respect to 
x over [-1,1] then using the orthogonality property (A.3) we obtain 

Gfcff( x jM*-3 Z 1 (1 - x 2 )* +2 c£ +2 \x) ($t?Wte = _<#+*>(*,.), 

7n-3Ajv-2 ( I — 

which enables us to write 


An-xTn-s 


w ‘ A«-jt.(*j)G5^J ,, (*j)G}Jil , (*>) 


, 2 < j < N - 2, 


(20) 


jince G^ A+2) ( X ) is a constant. Using the recurrence relation (A.6) with n - N 3 we write 
(17) in the form 

( 21 ) 


2 »+ 5 rnv + A-i)r(iv + A) 


I T ^ a/^V 1 / . 

w > = rN-3)!r(iv + 2A + ‘ir (1 -*?) 2 Gi}S ,, ( x j)G ( ^ ) ( x i) 


for 2 < j < N - 2. 

Representations for the boundary weights 
integrals (A.8) and (A.9). Consider u>i which, in view 
written in the form 


W\ and wn-i are found using the 
of (11), (12) (15) and (16), may be 


W\ 


1 

Gj}!?(-l)S((-» 



_ ,*)*(! _ *)*(! + 


( 22 ) 


6 


wrilepif in the form 0 + ^ ' ?!( ~ 1) ~ 4 ' The condithn ( A5 > to 


_ = HQiMMUi) /' 


Jji - x) J «(i + x) A +'G<J«>(x)rfx. 


4r(7V + A) 

The integral in (23) may be evaluated analytically using (A. 8) to give 

_ _ 2^+2p (A + 2 ) F (A + 3)^ _ 3) ; 

r~(7V + 2A + 2) • 


(23) 


(24) 


Similarly we can show that 

_ 2 2A+2 r(A + l )r(A + 3) (TV - 3)1 

^ (A + 3)r(7V + 2A + 2) {(A + 2)7V 2 + (A + 2)(2A — 1 )N — (4A 2 + 9A -f- 3)}, (25) 

and also 


w N—i — W\ , Wn~\ - —w t . 


(26) 


In the special case A - 0 the Gegenbauer polynomials coincide with the Legendre poly- 
nomials since ui 0 (x) = 1. The boundary weights are given by 


Wi = Wn_ 1 = 


Wi - -W N _i - 


8(2 TV 2 - 27V - 3) 

3(7V-2)(Ar- 1)7V(7V+ 1 ) ’ 
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and the interior weights by 


(TV — 2)(7V — 1)7V(7V + 1) ’ 


Wi — 


32(7V — 1)7V 


(27) 

(28) 

(29) 


' 2)(7V + + 2 ) 2 (T- xW&iUxj))* ’ 

f° r 2 < j < TV — 2. This form for the interior weights is derived using (A. 6) and (A. 10). 

When A = -1/2 the Gegenbauer polynomials are multiples of the Chebyshev polynomials 
T„(x) = cos(n cos **). In this case the boundary weights are given by 


W\ = Wjv_ 1 = 


_ 37r(3TV 2 - 6TV -+ 1) 


W\ = -WN - 1 = 

and the interior weights by 

tt(TV- 1)(TV - 2 ) 3 


10(TV — 2)(TV — 1)7V ’ 
37T 


4(TV — 2)(7V — \)N ’ 


Wj --- 


N (1 - 


, 2 < j < N - 2. 


(30) 

(31) 

(32) 


a„d H d aV T” S dOWn u an ex P, ression for » generalized psendospectral approximation (9) 

nd determined the weights in the corresponding quadrature rule we are now in a position 
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to write down the discrete problem corresponding to (4). The discrete variational problem 
corresponding to (4) is : Find tt;v £ Pn{~ 1>1) ^ ^a,o( ’ ) suc a 


a x (u N ,v N ) = {f,v N )\,d , V v N € P N {-^A) nH lo( 1 ’ 1 )’ 
where the bilinear form (., .) v is defined by 


(33) 


{f,g)x,d = N ^wj(xMx } )+ w ^f9)(xi)+(f9)(xN-i)]+M(f9y(xi)-(f9y(x N -i)}, (34) 

j= 2 

and Xj,2<j<N-2, are the zeros of GyJ” (*)• 

Theorem 3.1 The variational problem (33) is equivalent to the following collocation prob- 
lem: Finiu N € /M-1.1K' Hy-iA) such that 

rljf'to) = /(*j) , 2 < j < N - 2, (35) 

where xj, 2 < j < N - 2 are the zeros ofG^tl\x). 

Proof. The left-hand side, a*(u„,e„) is integrated by parts twice to give the following 
equivalent problem : Find un € flv(-U) n H A.o(~l. >) such that 


(u$\vn)\ = (/, V'ijxJ , 


(36) 


(38) 


for all v„ 6 ftr(-l. l) " H£o( 1 - 1 )' the , omia u h ( x ) , 2 < j < 

Now a basis for the space Pjv(-l,ljn«A,ot 1 ’^ arei11 F J (iw) ' , . , 

JV-2, defined in (10). These are used as test functions in (36). Since unUn y^lwehale 
HU-l, 1) and the quadrature rule (13) is exact for all polynomrals rn Pw-»(-l. >) 

(«<i”>, hi), = w.uWM ,2<j<N-2. (37) 

Further, from the definition (34), 

(/, hj)x,d = Wjf(xj) 2 < j < N - 2 , 
and therefore since > 0, for 2 < j < N - 2, we obtain (35) which completes the proof 

° f ‘we h“e“' analogous result to Theorem 2.1 for the discrete problem: 

Theorem 3.2 Let A satisfy -1 < A < 1. For any f € C°(-l,l), the problem (35) has a 
unique solution ujv € FW(— 1,1) H H%p( Ifl)- 

Bernard! and Maday (1991) establish the following error estimate: 

Theorem 3.3 Let A satisfy -1 < A < 1 . If the solution u of (29) belongs to ffj(- 1,1) for 
a real number a > 1, and if the data f is such that the function (1 ~*>f belongs to a space 
H >(- 1,1) for a real number p > 1/2, the following error estimate between the solutions of 

(29) and (38) is satisfied 

|| „ - uw |kr< C(N'-° || u |U +N || (1 - *¥/ |U). (39) 
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The collocation method (38) results in a system of equations for the values u, of uM 
^ e by° m ^ ’ 2 S ^ N ~ 2 ' The P seu dospectral collocation approximati™ is «L 

TV— 2 

u N (x) = Y, u jhj(x) , ( 40 } 

j = 2 ’ 

where h 0 (x) is given by (10) (cf. (9)). 

The generalization of the collocation method (35) for problems with inhomogeneous 

™”' tuchTat T 18 Stra ‘ ghtf ° rw " d - The nature <> f pseudospectral approximation 
( ) is such that inhomogeneous boundary conditions are satisfied exactly by simply insert- 
ing the specified values directly into (9). If H^- 1, 1) is the subspace of ^ TlVXch 
consis s o ose functions that satisfy the given inhomogeneous boundary conditions then 
we have the collocation problem : Find u N € P N (- 1, 1) n 1) such that 

u i N ) (x j ) = f( Xj ) , 2 < j < N — 2, 
where x j , 2 < j < N — 2 are the zeros of (?jvT^(x). 

4 Preconditioning 

Jq h uat C ions Cati ° n Pr ° b ' em (3S) ^ in the f ° rm of a linear ° { algebraic 

Au = b ( 41 ) 

where u is the vector of values of ti*(x) at the collocation points Xj , 2 < j < N - 2 b is the 

r^irir /(x) at these points and 4 is the in - 3)x <*-*>'■£* 

A J-hk-i = h ( k ,v) (xj), 2 < j,k < N - 2. 

he fourth order pseudospectral differentiation operator A has postive, real eigenvalues The 

ofr^Sr^lil Th Sh ° W n jn t TabIe 1 • I 1 " thiS Uble We 866 that the lar « est e 'S enva lue 
.. . bl e tbe sma llest eigenvalue is independent of N. Therefore since the 

.wlTof^r ° f 4 " 0(N ’ ] “ effiCiMt - necessary for^lcnl 

Orszag (1980) proposed a preconditioner for spectral methods based on a low-order finite 
mnH^ nCe appr ™ abon to the same differential operator. The advantages of such a pre- 

of the orid T h V *} 18 SParSG ’ eaSlly invertible and y ield s an inverse close to the inverse 
. S na s P ectral operator. Therefore we propose using a second-order finite difference 
ra or as our preconditioner. This requires the solution of a pentadiagonal system which 
m% be Performed very efficiently and stably using Gaussian elimination" " 
x . is ° r ~ 1 ~ 3 the second -order finite difference approximation to u(x) at the point 

u^(xi) k a iU (xi_ 2 ) + 6,u(*«-i) +' Ciu(xi) + div(x i+1 ) + e t u(x i+2 ) 

where 


a, = 


24 


( X *~ 2 x *)( x i- 2 - X,_ 1 )(x,_ 2 - X i+1 )(xi- 2 - X i+2 y 
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bi 

di 

Ci 


-24 


(x,_i - Xi)(Xi-2 - Zj_l)(Xj_l _ Z ‘+ 2 ) 

-24 

( Xi - X i+ i)(x;_2 - X i+ i)(x,_l - X, + l)(Xi+1 - Xi +2 ) 

24_ 

(Xj — Xi + 2)(Xi_2 — Xj+ 2 )(^t-l — Si+2)( x t+l — X, + 2 ) 
—(a* + bi + di + c t ). 


A similar formula h°^ s ^ 

(Ai\ The eigenvalues of H~'A are real and positive, ihe extreme eigen 

41 'i. ■ t uio 9 Again the smallest eigenvalue remains independent of the choice 

are shown in Table 2. Again 6 , „ r™ entr i es j n this table 

of N while the largest eigenvalue grows very si y • „ ldenw et a j (1984) 

(*%)<. We can see from Table 2 that they do indeed he between these boun . 

5 Analysis of the Multi-Domain Problem 

Given a fixed integer M we consider a partition of (-1.1) into M subintervals I m where 

J m = (d m , <i m +i), 

and the d m are M + 1 points in (-1,1) such that 

— 1 = d 0 < di < ■ • • < dM-i < djvf = I- 

supp ° rt on the interval /m ' 

We introduce the finite dimensional spaces 

Y n = {<f> € I 2 (-l,l) : ^ I /m£ fW(fm)}, 

where AT is some integer and PnW denotes the set of all 

equal to JV over A. In order to discretize the space H 0 {- 1, 1), let us imr 

Xn = Tat H 

Z N = {<t> € L 2 (— 1,1) '• |/ m € ^w(/m) H Hq ( Im , ) ) • 

The elements of are continuous and have continuous derivatives at the points «C, 1 < 
m < m - 1, and vanish along with their first derivatives at x - ±1. 


10 


concerned ThT T,?*" T'* ‘ he case A = 0 “ ^ as domain decomposition is 

concerned. This is the only value of A for which the weight function over (-1,1) is the 

same as e weight function over each of the subintervals I m 1 < m < M Throughout 

his section, in which A = 0 the zero subscript has been ddetedT For'exampl af ) i 

used synonymously with «( ) 0 and the corresponding norm notation has also been altered 

accordingly. We now define the discrete problem: Find u N 6 X N such that 


a(u N ,v N ) = Vv N 6 X N , 

where the bilinear form (., .) M is defined by 


M 


(fid)M = 53 ( f,g ) m , 


m= 1 


( 42 ) 


( 43 ) 


where 


N - 2 


g *>? f{x?)g{x?) + W?[(fg)(x?) + (/$)(:>#_,)] + W?[(f9)'(x?) - ( fg)’(x £_,)]. 

Lemma 5.1 For any real number a > 2 and for any <f> € H 2 {- 1, 1) n H*(- 1, 1) we have 

II <t> ~ *N fi 4 |fs < CN 2 ~° || $ || ff , 

where 7r 2 No is the orthogonal projection operator from H 2 {- 1, 1) onto P N (-1, l)fl tf 2 (-l, 1). 
Proof. See Bernardi and Maday (1991). □ 

Since /f 2 (— 1 1) is contained in C ,1 ([— 1, 1]) we can show that for any ^ e f/ 2 (— 1 1) there 
exists a cubic polynom.al «t„ such that * - & e »*(-!, l) a „ d for 4 teal number s > 0 

II * ||.< C II 4> 111 . 

S ^ M + *• fr0m *> W-h 1). So 

II ^ ^ IU = || (0 - - Jr* iO (0 - (f> 0 ) h 

< CN 2 ~ a \\ $- | )<r 

< CN*-' ||0|| r . 

We can easily verify that this operator satisfies 

(rU)(±l) = rf(±l),0r**)'(±l) = ^(±i). 

Theorem 5.1 There exists an operator P N from H 2 (- 1, 1) onto X N satisfying 


II 0 - || a < CN 2 ~ a || 0 ||„ 

for any function 0 € ff'(-l, 1) D tf 2 (-l, 1) * > 2. 


(44) 
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II ill' III 


Proof. We recall that for a general interval (a, b) there exists a projection operator from 
H 2 (a, b) onto P/v(a,6) satisfying 

|| w — |j H 2 (a,6) — CN 2 II w ||// <r (a,fc)» (45) 

n^w(a) = w(a), (t r^w)'(a) = w'(a), 

nNw(b) = w(b ), (7rAnuy(6) = u)'(6), (4^) 

fOF Let^ defile the projection operators T* m , for 1 < m < M, as being the projection 
operators from H 2 (I m ) onto P N {I m ). We deduce that the element defined on each I m 

y 7r 2 n 4>{x) = 7 r^ im ip(x), Vx € I m , 

is an element of Pjv(-1, 1) n // 0 2 (-l, 1) that satisfies due to (48) 

II if - lh< CN 2 ~° II if II. (48) 

Define J N f to be the Lagrange interpolating polynomial which interpolates the function 
/ at the N - 3 interior collocation points of the generalized Gauss quadrature rule on ( , )• 

Then Bernardi and Maday (1991) have shown that 

Lemma 5.2 For any real number p > 1/2 and for any <f> such that the function (\-x 2 )*<fe 
the following inequality holds 

|| (i _ x 2 )l(/ - J N f) || 0 < CN' /2 ~ p || (1 - x 2 ) 1 / IU • ( 49 ) 

Lemma 5.3 For any real number p m > 1/2 and for any f such that the function [(<*m+i - 
x)(x - d m )] 3/2 f 6 H pm (I m ), then 

z < CN^ 2 ~ Pm II [(«L+i - x)(x - d m )] 3/2 f ||iJptn(/ m ), 

v N tP N (iSnHZ(i m ) II VN (50) 

where (.,.)j m is the L 2 inner product on I m . 

Proof. The generalized Gauss quadrature rule on / m is exact for any polynomial in 
P 2 jv_ 3 (/ rn ) and so for any vn € Pv(4n ) H H 2 (I m ) we have 

{f,V N )l m -{f,VN)m = {f-JNf,V N )l m , 

where J N is the Lagrange interpolation operator at the N - 3 interior nodes of a generalized 
Gauss rule on the interval I m . We recall that the mapping w » w/[(d m+ 1 x){x d m ) j 

continuous from H 2 (I m ) into L 2 (I m ). Then we can write 

(f, VN ) Im - (f,v N ) m < C || [(cUi - *)(* - ^)l 3/2 (/ - J » f) ' 

Finally using Lemma 4.2 we obtain 

(| }U „) /m _ (f,V N ) m < CN ll2 ~ Pm II [(dm+l - *)(* - d m)f 2 f llHP»(/m)ll V » II ( 51 ) 

from which we deduce the result. O 


12 


Theorem 5.2 Let us suppose that the solution u to (32) belongs to H a ( — 1,1), for a real 
number a >2 and that [(rf m+1 - x)(x - d m )] 2 f £ H pm (I m ) for a real number p m > 1/2 for 
each m = 0, 1, . . . , M — 1. Then the following error estimate holds: 


A /-1 

II “ - “« h< C(N'-° II u II. + £ N II [ (<U, - *)(* - d m )] 3,2 f (52) 

m— 0 


Proof. Let us define u* — u - u° and u* N — un - u%, where u° and u% are piecewise cubic 
polynomials such that u* |/ m 6 H 2 (I m ) and u% | /m € P N (I m ) n H 2 (I m ), 0 < m < M - 1. 
ihen Proposition 2.1, together with (4) and (42), gives for any v N 6 Z N , 

II U N ~ V N H2 S C , ( a ( u /V — v Ni — Vn) 

= C(a(u - vn,u* n - vn) - (f,u* N -v N ) + (f,u* N - vn) m ), (53) 

from which we obtain 


u * ~ u *n II2 < C s inf 


k v n€%n 


+ sup (i' 1 ”") -(/■“■«)«) 


w N ez N 


W N 2 


■ 


We choose vjsj — and use Theorem 4.1 to show that 


(54) 


J”L II “* _ ” N II 2 - CN * ° II “* ll«> ^ > 2. (55) 

Since uiyv € Zn we have on each interval 7 m 

(/ » w N)l m ~ ( fi Wn)™. = (/ — J N f, w N)l m ■ 

We can also show that 


nm (f^N)~(f,W N ) M ^ 

&U P 71 n S sup — 

w N €Z N II ION \\2 m=0 ^6P w (/ m )nf/2 (/m ) || w N ||// 2 (/ m) 


(/» **>*)/„. - (/, WnU 


and therefore using Lemma 4.3 we may deduce that 

4“?„ S II [<*•« - *)(* - ^)] 3/ V l»-(W ■ (56) 

Since 


U- U N || 2 <|| U* - u* N || 2 4- II u° - || 2 

II u°-u° N h< C\\ u'-iil 


'^N II 2 , 


and 
then 

|| *i — Wyv || 2 < (7 || tt* - U* N || 2 . 

Finally using (54)-(56) we obtain the result. □ 

We now set up the collocation scheme for the domain decomposition problem. We define 
un £ Xn which interpolates data at the points x™ , 1 < j < N — 1, 1 < m < M by 

u N (x) = £ u Th?(x) + («')rET(x) + (u'«_,^_,(x) , x € /„, (57) 

i= 1 
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where 


(58) 




u N = u 


(u')iv = («')r +1 ,1 


Theorem 5.3 The variational problem (42) with the discrete inner product defined by (43) 
is equivalent to the following collocation problem : Find u/v € such that 

u%{x?) = /(xj 1 ), 2 < j < N — 2 , \ <m < M, (59) 

Un( x ? + 1 +)~ u n( x n- i-) = < +1 r(*? +1 +) + w %-A x N-i-)> 

1 < m < M - 1, ( 60 ) 

u'n( x ? +1 +) ~ u n( x n-i~) = -^r +lr « +1 +) “ w N-l r ( X N-l-) 

-l ff? +1 r\xA 1 +) -u%-i r '( x N-i-)i 
1 < m < M — 1, ( 61 ) 


and where r(x) = u$(x) — f{x). 

Proof. Let us examine a{u N , v N ) defined by (3). By linearity we may write the integral on 
the right-hand side of (3) as the sum of integrals over each subinterval I m for 1 < m _ 
Subsequent integration-by-parts twice gives 

fyf M — 1 

o(u W ,t>N) = E / U%(x )v N (x)dx - E {MMJM5-.) - K”«](*j5-l)}. < 62 ) 

Jim m = l 


where [/)(») s /(»+)-/(»-) denotes the jump at x = v in / 

We choose as our basis for the space X/v the polynomials h j (x), 2 < J S ^ - A 
1 < m < A/ - 1 and ^..(x), 1 < m < M - 1. The use of ‘ hese 

test functions in (42) with the discrete inner product given by (43) results in (59)-(61) whic 
completes the proof of the theorem. D 

Remark 1 Note that in view of the expressions for the weights given in (27) and (28), 


Wx = WN = 0(N~ 2 ) , m = -w N = 0(N 4 ) , as N oo , 


and therefore from (60) and (61) we can write 

«N(^ +1 +)-«^_ 1 -) = 0(iV- 2 ), 

<(< +1 +) -<(*»- !-) = 0(AT- 2 ) , 

as N -» oo. Thus we have second and third order continuity at the interface asymptotically, 
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6 The Biharmonic Problem in Two Dimensions 

Consider the biharmonic problem 


V^(x,y) = f(x,y), 

in ft , 

(63) 

H x ^v) = 9 \{x,y) , 

on r , 

(64) 

fa( x >y) = 92(x,y) , 

on r , 

(65) 


W er ^. , ~ [ > ) x an d T is the boundary of fl. Grisvard (1985) shows that 

proyided the boundary data satisfies certain compatibility conditions there exists ib h 6 H 2 ( SI) 
satisfying (64) and (65). Since we are primarily concerned with the domain decomposition 
pro em we only consider the case when the weight function is unity. The analysis for the 
single domain problem is thus greatly simplified. 

In order to write down the variational formulation of the problem (63)-(65) we define the 
bilinear form on H 2 (Cl) x H 2 { Q): ' 


«( 0 , <t>) = j jjy 2 xf>){V 2 <f>)dxdy . 


( 66 ) 

The biharmonic problem (63)-(65) is then equivalent to the following variational problem: 
Find if, 6 H 2 ( ft) such that (if, - i/, b ) € ff»(Q) an d 


4>) = (/, <t>) , for all <f> € H 2 (Q) , ( 67 ) 

where 

(/> 4>) = Jj^ f<j> dxdy . 

We see that 0 is a solution of the variational problem (67) if and only if X = 0 _ j,b [s a 
solution of the problem: Find xf> 6 H 2 ( ft) such that 


a(^, (j>) _ (/, <f>) - a(^ 6 , ^) } f or a u ^ € /f2(^) _ ( 6 g) 

eHiXon V6r ^m that , l 116 bili r ar form <-) defined b * (63) » continuous and 

J £ JJ- 2 (£|) 0 X an d hence that problem (71) has a unique solution in H£(fi,) for 

Let P N (n) be the space of algebraic polynomials of degree at most N in each co-ordinate 
direction, ihe collocation problem associated with (63)-(65) is- 
Find if> N € P N (tt) n H 2 (CL) such that 


V 4 Vw(z,J/) = f(x,y ) , 

(a:,y) € R, 

(69) 

^N{x,y) = gi(x,y ) , 

(x,y) eSUT, 

(70) 

fy>N ( v , . 

5n (^,y) -92{x,y) , 

(x, y)eSUT, 

(71) 


. (^,y)€T, 

(72) 
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where d/dn and d/dt represent differentiation normal and tangential to I\ respectively, the 
sets R, S and T are defined by 


S = {(^±l),(±l,6):2<*<N-2), 
T = {(±1,±1)}, 


. r (2) f/.'l = 0 2 < i < N — 2. There are a total of (JV + l) 2 linear equations for the 
Tn + untaow„;. The pension of P N (Sl) is (N + 1)*. The basis functions in 2 D are the 

tensor product of the one-dimensional basis functions given by (10) an ( )• 

We define the two-dimensional discrete inner product in an analogous way to (34) by 
applying the quadrature rule in each co-ordinate direction. So m the case when one of the 

functions or <j> belongs to //q( 0) we have 




N-2 N-2 

EE 


i=2 j = 2 


ui.uiji/ifC.O^fC.fj) • 


(73) 


Next we define the discrete bilinear form ajv(-> ■) by 

«jv(^j = (V 4 t/>,<^)jv ■ 


(74) 


Theorem 6.1 If there is a function € P„(n)n// 2 (ft) satisfying the boundary conditions 
(67) (68) then the collocation problem (69)-(7S) is eguiva lent to the vanaUonal problem. 
Find € Pv(ft) n ff 2 (ft) such that - V'fv) £ H o( n ) and 


a N (th ,,M = (/. Mn , /<”• » H f’n 6 fl*( n ) nfl 2( n ) 


(75) 


Proof. On each horizontal or vertical side of ft. and * are polynomials of degree 
AT satisfying N + 1 conditions and so are identical on T. The same argument applies to 

Lr normal derivatives and so (*„ - tfr) <E P«(ft) 0 W 0 J (ft). If we now choose *»{*, ») « 
" 2 < j,k < N - 2, then (75) implies (69)-(72). Conversely, smce these (IV - 3) 

polynomials form a basis for P„ 0 tfj(ft), (69)-(72) implies (75). □ 

Let us now turn our attention to the problem of domain decomposition, and for simplicity 
restrict onrsles for the moment to the case when ft is divided into two subdomains with 

7 ={(».»): -IS, SO- 

We define „ , . , t1 

ft 1 = {(x,y):-l<*<0, -1 < y < 1} , 

n a = {(x,y):0<*<l, -1 < y < 1} i 

and T k is the boundary of Q k for k = 1,2. Define the subspace V of H 2 ( fti) x H 2 {U 2 ) by 

, , dx}>' drP 2 , 

V = {* = (V’\V> 2 ) € H 2 (n a) x // 2 (D 2 ) : 0 = 4> , ^ ^ on ™ ’ 
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and the subspace Vq of V by 


^0 = {* = (V'W’ 2 ) € V : = 0 on r for m = 1,2 } . 

on ’ J 

We consider the bilinear form defined on VxVby 

a('P,$) = a 1 (V» 1 ,^ 1 ) + o 2 (V' 2 ,<A 2 ) , (76) 

where 

a k {tl> k ,<f> k ) = JJ q (VV‘)(VV)Wj, , 

We can show, using Green’s theorem, that if $ € V 0 then the above bilinear form may 
be written as 


«(*.*)= jj a (V'H>')<t>'ixdy + jj a (V'^)4?dxdy 

+ l ~ i ’ 2) % iy - 1 • (77) 

If there exists <P 6 € V satisfying (64)-(65) then the variational problem is: Find <P € V such 
that * - tf 6 e V 0 and ' 


«(*. *) = JJ ai ft'dxdy + JJ^ pcfdxdy , 
where f k is the restriction of / to 17*. 

The variational problem (78) is equivalent to the following interface problem: 


VV* = f, in ft*, fc=l, 2 , 

4> k = 9i , on rnr* , Jb = i, 2 , 

dip k 

= 92 , on rnr*, A: = 1 , 2 , 


0 1 = y> 2 , 


df_ 

dx 


Define the finite dimensional space VJv by 


5y> 2 

17 ’ on7 - 


v " = {* = (7\y> 2 ) e Pjv(fti) n tf 2 (ft x ) X P^(ft 2 ) n // 2 (ft 2 ) : 


and the subspace V)v,o of V)v by 


V- 1 = V> 2 , 


ay. 1 

5a: 


djP 

dx 


on 7} , 


(78) 


(79) 

(80) 

(81) 

(82) 


b)v,o {'P — (0 , y> ) € V^v : y> m = — — = 0 on T for m = 1, 2 } . 

In the case when one of or <P belong to V N , 0 we define a discrete inner product by 

(*>*)* = {4>\4> 1 Yn + (^ A 2 )n , 
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where 


N-2N-2 
»= 2 j= 2 

}= 2 

for Jk = 1,2, where d/dn is the normal derivative to ft*. The discrete bilinear form on 
Vn x Vn is defined by , . ,, 

a N (^,<I>) = aj v (^*) + 4(V’ ,<£ ) , 

where 

+ E w[^e>.6)^(°.6) - ^(0,f,)/(0,{,)] , 

i—2 


for * = 1,2. 

Theorem 6.2 If there is an element ^ b N € V N satisfying (80) and (81) then the variational 
problem: Find ^n E Vjv such that (f^N ^n) ^ 

a N (*N,*n) = (A*% + (/“■*% . (83) 

/or oil $ = (d 1 , d 2 ) £ V(v,o f 15 equivalent to the collocation problem: Find Wjv € that 

bl\ 1* II 


1 n( x ; 

k 


dfl)N_ 

dn 

d 2 ij) k 


N 


dndt 


dg 2 
at 


(x, y) e Rk , 

* = 1,2, 

(84) 

(x,y) e UTfc , 

* = 1,2, 

(85) 

(x,y) € S'*, U 7* , 

* = 1,2, 

(86) 

■ (x,y)€T*, 

* = 1,2, 

(87) 


|^*jjr - *ir)(0,&) = -«JU(WJr - /)(0,fc)J 

-w?[( v Vjv-/)(o,6-)] 

— W?[^(VV]v — /)(°»€i)] » 2 < j < N —2 , (88) 

C7X 

~(4>l - *'»)( o,fe) = -WfvVJ, - /)(»,&)] 

+ wj [(VVjv - /)(0, 6)1 , 2 < j < N - 2 , (89) 
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where 


R k = {({*,&) : 2 < i,j < N-2},k= 1 , 2 , 
= {U*,±l),(-l,6):2<*<7V-2}, 
s * = {(£?,±l),(l,&):2<t< N-2), 

T i = {(-1,±1),(0,±1)}, 

T 2 = {( 0 , ± 1 ), ( 1 , ± 1 )}, 


$ = (& - J )/ 2 - 6* = (1+ 60/2 , 2<i<N -2. 


Proof. We can show 
as our test functions the 


- 'f'/v) € V Nfi as in the proof of Theorem 5.1. If we now choose 
following: 


^{ x jy) — hk(2x -f 1 )hi(y) , <f> 2 (x,y) = 0 , 

<f>\x,y) = 0 , <j> 2 (x,y) = h k (2x-l)h,(y ), 

^ f x ’ y ) = h_N-\{2x + 1 )hi(y) , <f> 2 (x,y) = hi( 2x - l)h t (y) , 
• ^ = hN-i(2x + l)A/(j/) , <f> 2 (x,y) = A,(2x - , 


2 <k,l < N -2 , 
2 < k,l < N — 2 , 
2 < l < N — 2 , 

2 < 1 < N — 2 , 


(90) 


then we obtain immediately (84), (85), (88) and (89). Conversely, since these 2(N-3)(N- 2) 
test functions constitute a basis for V Nfit (84)-(89) implies (83). 


7 Numerical Results 


The quadrature rule (8) is used to compute approximations to the integrals 

(a) fl i u» l (x)e x cos(7rx)dx, 

0>) f-i w\(x)xe x dx, 


when A ~ 0 and A = -1/2. The errors in the quadrature rule are given in Tables 3 and 4 for 
integrals (a) and (b), respectively, for different' values of N. The quadrature rule evaluates 
the integrals accurate to machine precision for a value of TV as small as 17. 


7.1 1-D Problems 

Numerical solutions to the fourth-order model problem (1) are obtained when the exact 
solution is given by 


(a) u(x) = (1 - x 2 ) 2 sm( 7 rx), 

(b) u(x) = 1 -f sin(2irx). 


In example (a) the boundary conditions are homogeneous whereas for (b) we have inho- 
mogeneous boundary conditions. The differential equation is collocated at the generalized 
Legendre and Chebyshev nodes given by the zeros of (1 - x 2 )P£_ j(x) and (1 - z 2 )7^_j(x), 
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respectively. The error in the numerical solution is measured using weighted norms based 

on the corresponding generalized quadrature rule. The infinity norm r f Tabled 5 and 

the maximum pointwise error at the collocation points. These are displayed in Tables 5 and 

6 for examples (a) and (b), respectively, where we define 


N - 1 


1 2,u/ 


= [^^ J 2 +W 1 ((e' 1 ) 2 -(4_ 1 ) 2 ))] 

i= 1 


1/2 


,= max 
i <i<JV — i 


I e j I ’ 


and 


ej = u(xj) - Uff(xj) , 


= 1,2, iV — 1 , 


v. *L • - 1 < i < N — 1 are the generalized nodes. The usual exponential con- 

vergenceof^spectral "approximations' to^smooth solutions of differential equations is observed 

with accuracy to machine precision being obtained [when JV • For simplicity we 

T [-Ml and [ P 0,i] with 

common point * = 0. We solve again the model problems (a) and (b) using the collocation 
scheme (59)-(61). The corresponding error norms are shown in Tables 7 and 8, respe 
3? The mono-domain and two-domain spectral approximations converge exponentially 

defined in nonrectangular geometries in 2-D. 


7.2 2-D Problems 

Numerical solutions to the biharmonic equation are 
when the exact solution is given by 


obtained using the pseudospectral method 


(a) y) = (1 - x 2 )\\ - y 2 ) 2 sin{ny), 

(b) xl>{x,y) = (1 - * 2 ) 2 (1 - y 2 ) 2 sin{xx)sin(ry), 

(c) ip(x,y) = sin(2xx)sin(2xy). 


In examples (a) and (b) the boundary conditions are homogeneous whereas for (c) the 

is collocated at the Cartesian product of generalized Legendre nodes^ The w g d d 
infinity norms of the errors are shown in Table 9 for problems (b and (c). We see tha 
numerical solution correct to machine accuracy is obtained on a grid as coarse as 21 x 21. 
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In the case of domain decomposition we consider a partition of 0 into two rectangular 
subdomains [-1,0] x [—1,1] and [0, 1] x [-1,1] with common interface x = 0. Here we solve 
the model problem (a) using the collocation scheme (84)-(89). Both the mono-domain and 
two-domain pseudospectral approximations converge exponentially as expected. Again the 
two-domain approximation converges slower than the mono-domain approximation for the 
same total number of collocation points. This phenomenon was observed for 1-D problems 
too. In Figure 1 we give the contours of the approximation to the solution of problem (a) 
obtained using domain decomposition with N = 20. This figure is included to show the 
smoothness of the contours across the interface x = 0. 

7.3 2-D Problems in Nonrectangular Domains 

We extend the ideas developed in this paper to the solution of the Stokes problem for the 
How through an L-shaped channel. The flow geometry in this example is nonrectangular 
or which a standard single domain pseudospectral approximation is not applicable. The 
ability of pseudospectral methods to solve problems in this kind of geometry justifies the 
development of the theory of the multi-domain formulation considered earlier. The flow 
domain is divided into three rectangular subdomains as shown in Fig. 2. The stream 
function within each subdomain is approximated by a pseudospectral representation which 
interpolates values of the stream function at interior collocation points and values of the 
stream function and its normal derivative on the boundaries and subdomain interfaces. These 
representations are C l continuous across the subdomain interfaces. The unknowns in the 
pseudospectral approximations are determined from the collocation scheme derived from 
t e discrete variational formulation. This scheme results in C 3 continuous approximations 
asymptotically. 

If approximations of degree N are used in each direction in each subdomain then the 
collocation equations yield a system of (3 N - 5 )(N - 3) equations for the (3N -5)(N- 3) 
unknowns. A total of 2{N - 3) of these unknowns represent the values of the normal deriva- 
tives of at the interior nodes along the interfaces between subdomains and fi 2 and 
between subdomains 0 2 and H 3 . The remaining unknown values are the nodal values of xh 
at the interior and interface points of subregions Cl u and fi 3 . The collocation equations 
give rise to a linear algebraic system Au = b. The vector u contains the nodal values of tb 
and also the normal derivative of ^ at the interface nodes. The block tridiagonal structure 
ot the matrix A for the T-shaped domain is shown in Fig. 3. This system is solved using a 
Grout factorization subroutine from the NAG Library (1988). A more efficient direct solution 
echmque which takes account of the inherent matrix structure is the almost block diagonal 
solver of Brankin and Gladwell (1990) which has been used in spectral calculations by Kara- 

georghis and Phillips (1990,1991). However, this subroutine has not yet been incorporated 
into the present algorithm. 

The entry and exit lengths, a and b respectively, are chosen to be long enough to obtain 
AT y _ ^ ve flow * In Figs. 4 and 5 we show the contours of the stream function for 

_ ’ ~ 7 ’ f. = 1 Wlth a = -3 anci a = -5 > respectively. A small weak vortex is 

observed in the salient corner. Fully developed flow is reached within a channel width of the 
reentrant corner. 
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8 Conclusions 

Pseudospectral approximations to the solution of fourth-order elliptic partial differential 
equations are constructed using a collocation procedure based on the nodes of generalized 
Gaussian quadrature rules. Analytic expressions for the weights appearing in these qu 
ture rules are derived and their forms for the generalized Legendre and Chebyshev rules a 
given The equivalence between a discrete variational form of the differential problem wi 
suitably defined inner products and a collocation scheme is demonstrated when the coUo- 
cation points are chosen to be the zeros of certain ultraspherical polynomials. The natural 
choice of collocation points for fourth-order problems differs from the choice for second^ 
order problems, viz. the Gauss-Lobatto points. The usual convergence properties of spectral 

procedure based on the generalized Gauss-Legendre nodes is 
considered. Pseudospectral approximations which are automatically C - continuous at th 
subinterval interfaces are used to represent the solution An ^examination oi 
ing discrete variational problem results in an equivalent collocation method. The result g 
approximation is shown to be C 3 - continuous at the interfaces asymptotically, i.e. as the 
order of the approximations is increased in each subinterval. The scheme is analyzed and an 
error estimate is derived for the domain decomposed problem. 

For fourth-order problems in two dimensions we propose using a tensor product of th 
one-dimensional basis functions to represent the solution. The equivalence between the ^ 
cation method defined by collocating the differential equation on a grid formed by the tens 
product of the one-dimensional collocation points and a discrete variational formulatmn of 
the problem is described as well as the corresponding domain decomposition problem. It 
is intended to apply this collocation method to the solution of the Navier-Stokes equations 
in rectangularly decomposable domains using a stream function formulation even thoug 

simple variational principle does not exist for these equations. 

An application of this methodology to a biharmonic problem in a nonrectangular geom- 
etry is described. A single domain approach is.not feasible for this class of problem i unte 
one first transformed the original irregular domain to a simpler rectangular 
this would be cumbersome if it could be done at all since a transformation would need to be 

found for each new geometry. 
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TABLE 1 


Extreme eigenvalues of A 


N 


^mai 

ImajN* 

8 

0.3128+2 

0.4768+4 

0.2842-3 

12 

0.3128+2 

0.1091+6 

0.2537-3 

16 

0.3128+2 

0.1111+7 

0.2587-3 

20 

0.3128+2 

0.6788+7 

0.2652-3 

24 

0.3128+2 

0.2979+8 

0.2706-3 

28 

0.3128+2 

0.1039+9 

0.2750-3 

32 

0.3128+2 

0.3065+9 

0.2788-3 


TABLE 2 

Extreme eigenvalues of H~ l A 


N 

^mtn 

^mai 

8 

1.000 

3.312 

12 

1.000 

4.180 

16 

1.000. 

4.635 

20 

1.000 

4.915 

24 

1 1.000 

5.104 

28 

1.000 

5.241 

32 

1.000 

5.344 


TABLE 3 

Quadrature error in the approximation of w(x)e x cos(wx)dx 
for different weight functions 


N 

w(x) = 1 

w(x) = (1 — X 2 )" 1 / 2 

5 

0.497 -2 

0.689 -2 

9 

0.767 -10 

0.122 -8 

17 

0.300 -15 

0.710 -14 


TABLE 4 

Quadrature error in the approximation of f } , w(x)xe x dx 
for different weight functions 


N 

w(x) = 1 

w(x) = (1 — z 2 ) -1 /2 

5 

9 

17 

0.579 -5 
0.800 -14 
0.300 -15 

0.843 -5 
0.640 -14 
0.360 -14 
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TABLE 5 


Errors in the numerical solution of the model problem (1) 
with exact solution given by u(x) = (1 - x ) sm(irx) 


N 

A = 0 
II e II 2 . 11 / 

A = 0 

II e ll°° 

X = -1/2 
II elk. 

A = -1/2 

II e ll°° 

8 

1.387-2 

1.515-2 

3.041-2 

3.200-2 

12 

2.941-6 

2.954-6 

2.057-5 

2.783-5 

16 

3.077-10 

3.041-10 

6.355-9 

7.289-9 

24 

1.177-13 

1.329-13 

1.804-13 

2.018-13 


TABLE 6 

Errors in the numerical solution of the fourth-order problem 
with u(±l) = 1? du/dx(i: 1) — 2ir 
and exact solution u(x ) = 1 + sm(27r:r) 




A = 0 

A = 0 

A = -1/2 

A = -1/2 


N 

II e |1 

II e ll°o 

II e Ib.ui 

II e lloo 

. 

8 

0.236 

0.266 

0.450 

0.494 


12 

1.883-4 

1.845-4 

5.926-4 

7.976-4 


16 

1.181-7 

1.150-7 

8.848-7 

1.052-6 


24 

5.795-13 

6.701-7 

4.534-13 

5.190-13 


TABLE 7 


Errors in the 


numerical solution of the model problem (1) with exact solution given by 
= (1 - x 2 ) 2 sin(irx) using domain decomposition 


N 

A = 0 
II e II 2 , w 

A = 0 

II e ll°° 

A = -1/2 
II c \\2,w 

A = —1/2 

II e lloo 

6 

2.633-2 

2.024-2 

3.294-2 

2.747-2 

8 

1.490-3 

1.096-3 

1.890-3 

1.377-3 

12 

3.494-7 

2.642-7 

4.075-7 

3.073-7 

16 

1.293-11 

1.243-11 

1.254-11 

9.523-11 
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TABLE 8 


Errors in the numerical solution of the fourth-order problem 
with u(il) 1, du/dx{ il) 2it and exact Solution u(x) = 1 + sw(27rx) using domain 

decomposition 


N 

A = 0 

II e Ik* 

A = 0 
II e |L 

A = -1/2 

II « Ik* 

A = -1/2 

II e |L 

6 

0.326 

0.230 

0.476 

0.331 

8 

2.794-2 

2.080-2 

3.580-2 

2.632-2 

12 

3.221-5 

2.430-5 

3.757-5 

2.839-5 

16 

7.491-9 

5.677-9 

7.961-9 

6.058-9 


TABLE 9 

Errors in the numerical solution of the biharmonic problem (63)-(65) 
with exact solutions given by (b) and (c) 



Problem (b) 

; Problem (c) 

N 

1 ^ 2,u> 

e oo 

1 e 2,* 

I e I oo 

6 

0.347 

0.512 

1.467 

2.483 

8 

5.145-3 

3.779-2 

5.831-2 

8.238-2 

12 

4.667-6 

9.413-5 

1.879-4 

4.925-4 

16 

5.631-10 

7.836-9 

3.856-8 

7.437-8 

20 

2.935-13 

4.291-12 

9.648-11 

1.502-10 
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Figure 1. Contour plots of tl>(x,y) for problem (a) when N -20 , using domain 
decomposition and the generalized Legendre pseudospectral method. 
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Figure 4. Contours of 1>(z,y) for a = -3, b = 7,c = 1 and N = 16. 
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A Properties of Ultraspherical Polynomials 

The ultraspherical or Gegenbauer Polynomials are the solutions of the differential equation 

[wa+i(*)G 4 a >'(*)] # + n(n + 2A + l)u> A (x)C?( A >(x) = 0, (1) 

that are bounded at x = -h 1 } where 


w x (x) = (1 — ar 2 ) A , A > -1. 

They are orthogonal with respect to w x (x) over the interval [-1,1] : 

]_ x M^)G^\x)G^\x)dx = 7m 6 m , n , 

where 

2 2 ^[r(m + A + l )l 2 
(2m + 2A + l)m!r(m + 2A + 1) ’ 

and T is the gamma function. At x = ±1 , G< A >(x) satisfies the condition 

G* A) (±1) = ( ± l)n£(!i±A±i) 
n v ' k } n!r(A + l) ' 

The ultraspherical polynomials may be generated using the recurrence relation 
(n + l)(n + 2A + l)Gi A + ) 1 

= (2n + 2A + l)(n + A + 1)*G<, A > - (n + A)(n + A + 

Gr(x) = 1, G'i A) (x) = ( A + l)x. 

The leading coefficient, A n , of G^ A) (x) is given by 


( 2 ) 

( 3 ) 

( 4 ) 

( 5 ) 


( 6 ) 


A _ 1 T(2n + 2A + 1) 
2"n!r(n + 2A -+ 1) 


( 7 ) 


We have the following integrals involving ultraspherical polynomials (Erdelyi (1954), p.284) 

r (1 - s) A (l 4 - xYG^\x)dx = 2A+<7+lr ( g + l)r(A + n + l)r(q - A + 1) 

J n W n!T(<T — A — n + 1)T(A + cr + n + 2) ’ (8) 

( 9 ) 


/ (1 — xf{ 1 + xfG^\x)dx = - ^ + l)r(A + n + l)r(A - + n) 

J ~ l n!r(A - a)T(\ + cr + n + 2) 

where A, a > — 1. 

The ultraspherical polynomials satisfy the recursion relation 


(! - x 2 )Gi A )'(x) = ~nxG^\x) + (n + A)^, (x). 


( 10 ) 
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